The proportion of time an animal is in a feeding behavioral state.
Process Model
\[Y_{i,t+1} \sim Multivariate Normal(d_{i,t},σ)\]
\[d_{i,t}= Y_{i,t} + γ_{s_{i,g,t}}*T_{i,g,t}*( Y_{i,g,t}- Y_{i,g,t-1} )\]
\[ \begin{matrix} \alpha_{i,1,1} & \beta_{i,1,1} & 1-(\alpha_{i,1,1} + \beta_{i,1,1}) \\ \alpha_{i,2,1} & \beta_{i,2,1} & 1-(\alpha_{i,2,1} + \beta_{i,2,1}) \\ \alpha_{i,3,1} & \beta_{i,3,1} & 1-(\alpha_{i,3,1} + \beta_{i,3,1}) \\ \end{matrix} \] \[logit(\phi_{Behavior}) = \alpha_{Behavior_{t-1}} \] The behavior at time t of individual i on track g is a discrete draw. \[S_{i,g,t} \sim Cat(\phi_{traveling},\phi_{foraging},\phi_{resting})\]
Dive information is a mixture model based on behavior (S)
\(\text{Average dive depth}(\psi)\) \[ DiveDepth \sim Normal(dive_{\mu_S},dive_{\tau_S})\] \[ DiveDuration \sim Normal(duration_{\mu_S},duration_{\tau_S})\]
Dive Profiles with Argos timestamps
## # A tibble: 11 x 4
## Animal n argos dive
## <int> <int> <int> <int>
## 1 131111 548 173 375
## 2 131115 1029 179 850
## 3 131116 1932 457 1475
## 4 131127 9275 2495 6780
## 5 131128 112 52 60
## 6 131130 1165 151 1014
## 7 131132 2292 589 1703
## 8 131133 8497 2078 6419
## 9 131134 2098 653 1445
## 10 131136 6844 1970 4874
## 11 154187 1866 486 1380
## [1] 35658
## # A tibble: 11 x 2
## Animal n
## <int> <int>
## 1 131111 454
## 2 131115 1020
## 3 131116 1918
## 4 131127 9047
## 5 131128 106
## 6 131130 511
## 7 131132 1713
## 8 131133 8193
## 9 131134 2018
## 10 131136 6721
## 11 154187 1847
## # A tibble: 11 x 2
## Animal Tracks
## <int> <int>
## 1 131111 1
## 2 131115 1
## 3 131116 2
## 4 131127 11
## 5 131128 1
## 6 131130 1
## 7 131132 1
## 8 131133 4
## 9 131134 1
## 10 131136 3
## 11 154187 1
sink(“Bayesian/ThreeState.jags”) cat(" model{
#from jonsen 2016
pi <- 3.141592653589
#for each if 6 argos class observation error
for(x in 1:6){
##argos observation error##
argos_prec[x,1:2,1:2] <- argos_cov[x,,]
#Constructing the covariance matrix
argos_cov[x,1,1] <- argos_sigma[x]
argos_cov[x,1,2] <- 0
argos_cov[x,2,1] <- 0
argos_cov[x,2,2] <- argos_alpha[x]
}
for(i in 1:ind){
for(g in 1:tracks[i]){
## Priors for first true location
#for lat long
y[i,g,1,1:2] ~ dmnorm(argos[i,g,1,1,1:2],argos_prec[1,1:2,1:2])
#First movement - random walk.
y[i,g,2,1:2] ~ dmnorm(y[i,g,1,1:2],iSigma)
###First Behavioral State###
state[i,g,1] ~ dcat(firstmove[]) ## assign state for first obs
#Process Model for movement
for(t in 2:(steps[i,g]-1)){
#Behavioral State at time T
phi[i,g,t,1] <- Traveling[state[i,g,t-1]]
logit(phi[i,g,t,3]) <- alpha_hours[state[i,g,t-1]] + beta_hours[state[i,g,t-1]]* cos((2*pi*hours[i,g,t])/(24)) + beta2_hours[state[i,g,t-1]] * sin((2*pi*hours[i,g,t])/24)^2
phi[i,g,t,2] <- 1-(phi[i,g,t,1] + phi[i,g,t,3])
state[i,g,t] ~ dcat(phi[i,g,t,])
#Turning covariate
#Transition Matrix for turning angles
T[i,g,t,1,1] <- cos(theta[state[i,g,t]])
T[i,g,t,1,2] <- (-sin(theta[state[i,g,t]]))
T[i,g,t,2,1] <- sin(theta[state[i,g,t]])
T[i,g,t,2,2] <- cos(theta[state[i,g,t]])
#Correlation in movement change
d[i,g,t,1:2] <- y[i,g,t,] + gamma[state[i,g,t]] * T[i,g,t,,] %*% (y[i,g,t,1:2] - y[i,g,t-1,1:2])
#Gaussian Displacement in location
y[i,g,t+1,1:2] ~ dmnorm(d[i,g,t,1:2],iSigma)
}
#Final behavior state
phi[i,g,steps[i,g],1] <- Traveling[state[i,g,steps[i,g]-1]]
logit(phi[i,g,steps[i,g],3]) <- alpha_hours[state[i,g,steps[i,g]-1]] + beta_hours[state[i,g,steps[i,g]-1]]* cos((2*pi*hours[i,g,steps[i,g]])/(24)) + beta2_hours[state[i,g,steps[i,g]-1]] * sin((2*pi*hours[i,g,steps[i,g]])/24)^2
phi[i,g,steps[i,g],2] <- 1-(phi[i,g,steps[i,g],1] + phi[i,g,steps[i,g],3])
state[i,g,steps[i,g]] ~ dcat(phi[i,g,steps[i,g],])
## Measurement equation - irregular observations
# loops over regular time intervals (t)
for(t in 2:steps[i,g]){
# loops over observed locations within interval t
for(u in 1:idx[i,g,t]){
zhat[i,g,t,u,1:2] <- (1-j[i,g,t,u]) * y[i,g,t-1,1:2] + j[i,g,t,u] * y[i,g,t,1:2]
#for each lat and long
#observed position
argos[i,g,t,u,1:2] ~ dmnorm(zhat[i,g,t,u,1:2],argos_prec[argos_class[i,g,t,u],1:2,1:2])
#for each dive depth
#dive depth at time t
dive[i,g,t,u] ~ dnorm(depth_mu[state[i,g,t]],depth_tau[state[i,g,t]])T(0,)
#Assess Model Fit
#Fit dive discrepancy statistics
eval[i,g,t,u] ~ dnorm(depth_mu[state[i,g,t]],depth_tau[state[i,g,t]])T(0,)
E[i,g,t,u]<-pow((dive[i,g,t,u]-eval[i,g,t,u]),2)/(eval[i,g,t,u])
dive_new[i,g,t,u] ~ dnorm(depth_mu[state[i,g,t]],depth_tau[state[i,g,t]])T(0,)
Enew[i,g,t,u]<-pow((dive_new[i,g,t,u]-eval[i,g,t,u]),2)/(eval[i,g,t,u])
}
}
}
}
###Priors###
#Process Variance
iSigma ~ dwish(R,2)
Sigma <- inverse(iSigma)
##Mean Angle
tmp[1] ~ dbeta(10, 10)
tmp[2] ~ dbeta(10, 10)
tmp[3] ~ dbeta(10, 10)
# prior for theta in 'traveling state'
theta[1] <- (2 * tmp[1] - 1) * pi
# prior for theta in 'foraging state'
theta[2] <- (tmp[2] * pi * 2)
theta[3] <- (tmp[3] * pi * 2)
##Move persistance
# prior for gamma (autocorrelation parameter)
##Behavioral States
#gamma[1] ~ dbeta(3,2)
#dev ~ dunif(0,0.5) ## a random deviate to ensure that gamma[1] > gamma[2]
#gamma[2] <- gamma[1] * dev ## 2d movement for foraging state
#dev2 ~ dunif(0,0.5) ## a random deviate to ensure that gamma[1] > gamma[3]
#gamma[3] <- gamma[1] * dev2 ## 2d movement for resting state
gamma[1]<-0.7
gamma[2]<-0.3
gamma[3]<-0.3
#Temporal autocorrelation in behavior - remain in current state
Traveling[1] ~ dbeta(1,1)
Traveling[2] ~ dbeta(1,1)
Traveling[3] ~ dbeta(1,1)
#Temporal autocorrelation in behavior - transition Foraging
alpha_hours[1] ~ dnorm(0,0.386)
alpha_hours[2] ~ dnorm(0,0.386)
alpha_hours[3] ~ dnorm(0,0.386)
#Effect of time of day on transitioning to resting
beta_hours[1] <- 0
beta_hours[2] ~ dnorm(0,0.386)
beta_hours[3] ~ dnorm(0,0.386)
#Transition from travel to resting
beta2_hours[1] <- 0
beta2_hours[2] ~ dnorm(0,0.386)
beta2_hours[3] ~ dnorm(0,0.386)
#Probability of initial behavior
firstmove ~ ddirch(rep(1,3))
#Foraging dives are deepest
depth_mu[2] <- 0.150
depth_mu[1] <- 0.02
depth_mu[3] <- 0.015
#depth and duration variance
depth_tau[1] <-2000
depth_tau[2] <- 300
depth_tau[3] <- 5000
##Observation Model##
##Argos priors##
#longitudinal argos precision, from Jonsen 2005, 2016, represented as precision not sd
#by argos class
argos_sigma[1] <- 11.9016
argos_sigma[2] <- 10.2775
argos_sigma[3] <- 1.228984
argos_sigma[4] <- 2.162593
argos_sigma[5] <- 3.885832
argos_sigma[6] <- 0.0565539
#latitidunal argos precision, from Jonsen 2005, 2016
argos_alpha[1] <- 67.12537
argos_alpha[2] <- 14.73474
argos_alpha[3] <- 4.718973
argos_alpha[4] <- 0.3872023
argos_alpha[5] <- 3.836444
argos_alpha[6] <- 0.1081156
}"
,fill=TRUE)
sink()
## user system elapsed
## 2464.286 571.552 9275.116
Reconstructing the transition matrix
## user system elapsed
## 0.573 0.131 459.074
Where does the 3d predict foraging that the 2d misses?
Where does the 3d predict resting that the 2d misses?
Where does the 2d predict foraging that the 3d refines to traveling?
Where does the 2d predict foraging that the 3d refines to resting?
The goodness of fit is a measured as chi-squared. The expected value is compared to the observed value of the actual data. In addition, a replicate dataset is generated from the posterior predicted intensity. Better fitting models will have lower discrepancy values and be Better fitting models are smaller values and closer to the 1:1 line. A perfect model would be 0 discrepancy. This is unrealsitic given the stochasticity in the sampling processes. Rather, its better to focus on relative discrepancy. In addition, a model with 0 discrepancy would likely be seriously overfit and have little to no predictive power.
## # A tibble: 1 x 2
## `mean(E)` `var(Enew)`
## <dbl> <dbl>
## 1 0.1259669 0.008264216